Energy Landscape of the Finite-Size Mean-field 3-Spin Spherical Model 
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We study the 3-spin spherical model with mean-field interactions and Gaussian random couplings. 
For moderate system sizes of up to 20 spins, we obtain all stationary points of the energy landscape by 
means of the numerical polynomial homotopy continuation method. On the basis of these stationary 
points, we analyze the complexity and other quantities related to the glass transition of the model 
and compare these finite-system quantities to their exact counterparts in the thermodynamic limit. 
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I. INTRODUCTION 

In this work we analyze key characteristics of the en- 
ergy landscape of the spherical p-spin glass model (pSM) 
[1, 2] by means of the numerical homotopy continuation 
method. This method is capable of finding all stationary 
points of the energy function. However, since the number 
of stationary points grows exponentially with the system 
size TV, one is restricted to rather small systems. The in- 
terest in the pSM resides in its connection with the theory 
of the glass transition. Although it is, strictly speaking, 
a spin glass model with disordered quenched couplings, 
its physics resembles in many respects that of structural 
glass formers and has been a key model in the study of 
the glass transition (for a review see [3]). 

Its landscape properties have been extensively stud- 
ied in the large-iV limit where, owing to its mean-field 
character, the model can be solved exactly. The number 
of stationary points of the so-called Thouless- Anderson- 
Palmer (TAP) free energy has been computed and char- 
acterized in Refs. [2, 4-6]. A key quantity is the com- 
plexity 



S(e) 



N 



InAA(e) 



(1) 



as a function of the energy density e, where Af is the 
number of stationary points of the TAP free energy. Be- 
sides the number of stationary points, also their index 
/ turns out to be relevant, i.e., the number of negative 
eigenvalues of the Hesse matrix evaluated at a stationary 
point. For this reason it is convenient to also introduce 
the complexities 



EKe) = -InA/}(e) 



(2) 



where Afj is the number of stationary points of the TAP 
free energy with index /. The picture that emerges is that 
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the pSM has a high-energy regime characterized by the 
presence of exponentially many stationary points with 
large index. Asymptotically in the large- limit, the 
complexities E/ are equal for all / = 0, . . . , A^ — 1. This 
regime corresponds to a high-temperature, or "simple liq- 
uid" , regime. Below a certain threshold energy eth the 
number of minima becomes exponentially dominant over 
higher-index stationary points, and the system's dynam- 
ical behavior is governed by states near the bottom of the 
landscape [7] . Interestingly, eth coincides with the energy 
at which a dynamic singularity is present, in the sense 
that relaxation times diverge in the large- A^ limit. This 
dynamic singularity corresponds to the so-called mode 
coupling transition and is present in many glass form- 
ing models [8-10]. In the minima-dominated regime, for 
large A^, the complexity of minima gives a finite contri- 
bution to the free energy of the system, called the con- 
figurational entropy. 

In mean-field models like the one we are considering in 
the present article, divergent energy barriers lead to true 
ergodicity breaking at eth- Although this is known not 
to be true anymore in more realistic models with short- 
range interactions, several short-range models still show 
signatures of a change of regimes in the energy landscape 
near the mode coupling transition. One such signature is 
the vanishing of the mean index density i near the mode 
coupling transition, where i(e) is obtained by averaging 
the index density i = I/N of the Hesse matrix over all 
stationary points with energy e [9, 10]. In real systems it 
is expected that relaxation in the low energy regime will 
be driven by activation over barriers, due to the domi- 
nance of minima in the energy landscape probed by the 
system's dynamics. The complexity of minima stays pos- 
itive in a restricted energy window below the threshold, 
diminishing continuously until a second characteristic en- 
ergy below which the complexity changes sign. In the 
mean field picture this means that, for energy densities 
e < ec, the number of minima becomes exponentially 
small (instead of exponentially large) for large N. In- 
terestingly, this energy coincides with the point at which 
a replica calculation gives a replica symmetry breaking 
phase transition in the thermodynamic limit. At this 
point a true phase transition, hallmarked by a single 



2 



step of replica symmetry breaking, occurs and the model 
enters a spin glass phase (or "ideal glass phase" in the 
structural glass terminology [4, 11]). 

From a landscape perspective on thermodynamic prop- 
erties of the model, the pSM has the striking feature that 
the equations defining the stationary points of the TAP 
free energy (which is a mean-field free energy function of 
local magnetizations) are formally identical to the equa- 
tions determining the stationary points of the Hamilto- 
nian [2, 5]. Therefore, in the thermodynamic limit, many 
temperature-dependent properties of the model can be 
inferred directly from an analysis of the energy land- 
scape. In general (i.e., for models other than the pSM) 
such a duality of TAP free energy stationary points and 
Hamiltonian stationary points does not hold. The ex- 
trapolation of energy landscape characteristics (based on 
the Hamiltonian) to the thermal behavior of glassy sys- 
tems (governed by the free energy) must be done with 
considerable care. 

The energy landscape approach has greatly improved 
the understanding of the glass transition, and much of 
this success is due to the exact results for landscape- 
related quantities of the pSM in the thermodynamic 
limit. Nonetheless, many questions related to the way 
that barriers and timescales scale with size are relevant 
for real systems and are still poorly understood. The nu- 
merical homotopy continuation method, by virtue of its 
capacity of determining all stationary points of a complex 
polynomial, opens up a way of obtaining exact finite- 
system properties of the energy landscape of the pSM for 
system sizes up to = 20 (i.e., far from the thermody- 
namic limit). The main goals of the present study are (1) 
to characterize the landscape of a small cluster of spins, 
which can be considered as a prototype of a real, large- A'^ 
glass model, and approach in a systematic way the emer- 
gence of complexity as known in the large-A^ case, and 
(2) to test how far the validity of the thermodynamic- 
limit results for the pSM extends to finite and even small 
system sizes. 



II. THE MEAN-FIELD p-SPIN SPHERICAL 
MODEL 

The p-spin spherical model is defined by the Hamilto- 
nian 

N 



^ ii,Z2,...,ip — 1 

subject to the spherical constraint 



N 



(4) 



that restricts the configurations a = (cti, . . . , ctat) e 

to an (A^— l)-spherc of radius ^/N (hence the name of the 

model). The configuration space of the pSM is therefore 



effectively (A'' — l)-dimensional. The coupling constants 
Jii,...,ip are Gaussian random variables with zero mean 



and standard deviation ^/p\J2NP^. The Hamiltonian 
(3) consists of p-spin interactions between all possible 
groupings of p lattice sites. For such a fully-connected 
model the mean-field approximation is known to be ex- 
act, and for this reason the model itself is called the 
mean-field pSM. At least qualitatively, the phenomenol- 
ogy of the pSM is rather insensitive to the value of the 
integer parameter p, as long as p > 3. For this reason we 
will restrict our numerical study to the case p = 3 where 
the Hamiltonian of the 3-spin spherical model (3SM) is 
given by 



N 



i.j,k—l 



(5) 



In order to compute stationary points of the Hamilto- 
nian (3) restricted to the spherical configuration space by 
the constraint (4), it is convenient to resort to the formal- 
ism of Lagrange multipliers. We define the Lagrangian 



L 




N 



(6) 



where ao is a Lagrange multiplier. The stationary points 
of the constrained system are then the solutions of the 
system of equations 



dL 
da. 



0, 



i = 0,...,N. 



For the 3SM these equations are 



N 

Ji,j,kCrj(Jk + 6eCTj = 0, 



(7) 



(8) 



where e = H/N is the energy density at a stationary 
point [5]. The stability of a stationary point is governed 
by the Hesse matrix 



T~{-lm — 



daid(Tr, 



N 



(9) 



k=l 



with l,m = 1,...,A^, where denotes Kronecker's 
symbol. The matrix "H as defined in (9) is an A?^ x A^ ma- 
trix. Evidently, this cannot be the correct Hesse matrix 
of the 3SM: as mentioned earlier, the spherical constraint 
(4) restricts the configuration space of the model to an 
(A^ — l)-dimensional sphere. Hence the correct Hessian 
H, i.e., the one constrained to the spherical configura- 
tion space, should be an (A^ — 1) x (A^ — 1) matrix. This 
can be accounted for by eliminating from T-l the eigendi- 
rection perpendicular to the surface of the configuration 
space manifold. For a spherical configuration space such 
a perpendicular eigenvector is radial for any given sta- 
tionary point. It is possible to prove that the eigenvalue 



3 



corresponding to this radial eigenvector is equal to Be for 
the 3SM [6]. Knowledge of this radial eigenvalue allows 
us to easily obtain the index of the correct (constrained) 
Hessian % from the unconstrained one, 



/(H(a^)) = /(H(a^))-e(-e), 



(10) 



where = {<7\, . . . , a^^) is a solution of (8) and e = 
H{a^)/N is the corresponding energy density. Similarly, 
the determinant of the constrained Hessian can be ob- 
tained from "H, 



detCH) 



detCH) 
3e 



for e^O. 



(11) 



III. 



NUMERICAL POLYNOMIAL HOMOTOPY 
CONTINUATION METHOD 



The stationary points of the 3SM are the solutions of 
the A'^ coupled nonlinear equations in (8). Since an ana- 
lytic solution of these equations is unlikely to be feasible, 
we resort to a numerical technique. Due to the polyno- 
mial character of (8), the so-called numerical polynomial 
homotopy continuation (NPHC) method [12] is particu- 
larly suitable. This method has the virtue of being able 
to find all the solutions of a given system of polynomial 
equations, and has been applied in the past to a variety 
of problems in particle theory and statistical mechanics 
[13-24]. A drawback of the NPHC method is that it is 
restricted to fairly small system sizes. A detailed descrip- 
tion of the method can be found in [15] or in Sec. HI of 
[19]. 

Running the NPHC method on a computer, we con- 
sider a configuration (cti, . . . , ctat) to be a numerical solu- 
tion of (8) if it satisfies these equations with tolerance not 
exceeding 10~^°. The NPHC method searches solutions 
in and; since we are interested only in real solutions, 
we have to dispose of all the complex ones as a last step. 
We consider a solution to be real if the absolute value of 
the imaginary part of each component ai does not exceed 
the tolerance 10~^. We have checked that the choice of 
this tolerance does not affect the number of real solutions 
found, i.e., this number is robust for the problem at hand. 
All stationary points can be further refined to arbitrary 
numerical precision if necessary. Moreover, we have used 
an efficient implementation of Smale's a-theory [25] to 
certify that each of the numerical solutions is indeed in 
the quadratic convergence of the associated actual solu- 
tion (see [26] for an introduction of the method from the 
potential energy landscape point of view). 



IV. RESULTS 

By means of the NPHC method we computed all sta- 
tionary points of the 3SM for system sizes N between 14 
and 20 and, for each value of N, ten different realizations 
of the random couplings were used. In the Introduction 
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FIG. 1. (color online) Top: Disorder average E(A/') of the 
total number A/" of complex (red squares) and real (blue dots) 
stationary points on a logarithmic scale, plotted vs. the system 
size A^. Both numbers grow exponentially with A^. Bottom: 
Logarithmic density of E(A/') for the real stationary points, 
plotted vs. the inverse system size (blue dots). The analytic 
result in the thermodynamic limit is indicated by a red square. 



we had outlined the relevance of the stationary points 
and their indices for the relaxational properties of glassy 
systems. In the following sections we will analyze for 
the finite-size 3SM the most relevant properties, includ- 
ing complexities, mean indices, and the Hesse determi- 
nant, and compare them to analytic predictions for their 
infinite-system counterparts. 



A. Number of stationary points 

In Fig. 1 (top) wc show the total numbers of complex 
and real stationary points, respectively, averaged over the 
disorder samples. The linear increase in this logarithmic 
plot nicely illustrates that the total number of stationary 
points grows exponentially with N. While this is true in 
many models, the relevance of this behavior in the 3SM 
model is that not only the energy function, but the free 
energy itself has the same structure of stationary points. 
This is in stark contrast with the majority of statisti- 
cal systems for which the number of stationary points of 
the free energy in order parameter space is usually much 
smaller. This is intimately connected to the glassy be- 
havior of the 3SM at low temperatures/energies. The 
large number of complex stationary points is the limiting 
factor for applying the NPHC method to larger system 
sizes. 

In the large- limit, an analytic expression is known 



4 



for the disorder average E(A^) of the total number of 
real solutions, 



lim llnE(AA) = iln2, 

N^oo iV Z 



(12) 



(Eq. (2.20) of Ref. [27]). In Fig. 1 (bottom) we plot 
hi'E{Af)/N vs. 1/N for our finite-system data. Although 
not an exact match, the numerical data tend nicely to- 
wards the exact limiting value, even for the very small 
system sizes considered. 



B. Complexity vs. energy 

Resolving the number of stationary points with respect 
to the energy, we obtain the complexity S(e) as intro- 
duced in (1). Although Eq. (1) captures the essential 
idea of the complexity, a technical difficulty is hidden 
in the number Af{e) of stationary points at energy den- 
sity e: Rigorously speaking, Af{e) will be zero almost ev- 
erywhere. To obtain a well-defined quantity, some kind 
of binning would be necessary, for example by defining 
A/'(e) to be the number of stationary points in the inter- 
val [e, e -I- A] for some small A > 0. A drawback of this 
approach is that an arbitrary parameter A is introduced, 
and the effect of the choice of A is difficult to control. 

This problem can be circumvented by resorting to the 
cumulative complexity 



where 

M{e) 



r(e) = llnA^(e), 



(13) 



(14) 



is the number of stationary points with energies not 
larger than e. The cumulative complexity r(e), though 
discontinuous, is a well-defined, monotonically increasing 
function for all finite system sizes N. Also in analyti- 
cal calculations the cumulative complexity T is usually 
preferred to its non-cumulative counterpart S, and ex- 
act analytic expressions are known for the pSM in the 
thermodynamic limit (Eq. (2.15) of [27]). The cumula- 
tive complexity is a self-averaging quantity in the large- A'^ 
limit, but since our system sizes are far from that limit, 
we calculated quenched averages 



(r)c 



1 

N 



E{lnM) 



as well as annealed averages 



1 



(r)a = -lnE(X), 



(15) 



(16) 



where E denotes averaging over 10 disorder realizations 
[5]. In Fig. 2 we show (r)a for various system sizes N, to- 
gether with the thermodynamic limit result. Although, 




FIG. 2. (color online) The annealed average of the cumula- 
tive complexity F plotted vs. the energy density e for various 
system sizes A'^. For comparison, the analytic large-A'^ result 
is shown as a black line. 



for the small N considered, the finite-system data are 
not close to the thermodynamic limit result for most e, 
the trend we observe is correct: For energies e less than 
about —0.55 the finite-system data lie below the infinite 
system result, but become larger as N grows; for e larger 
than about —0.45 the finite-system data lie above the in- 
finite system result, but decrease with increasing N. The 
quenched average (F)q (not shown) turns out to be very 
similar to (F)a and the two are almost indistinguishable 
on the scale of Fig. 2 for most e. The only visible differ- 
ence is that, for energies smaller than the ground state 
of any one of the disorder realizations, (F)q is ill-defined 
(or — oo). As a consequence, the quenched average (F)q 
is nonnegative, and convergence to the negative values 
of the analytic result in the thermodynamic limit cannot 
occur. For this reason the annealed average (F)a appears 
more suitable for such a comparison. 

In particular, the critical energy Cc, i.e. the energy 
at which the complexity goes to zero in the thermody- 
namic limit, can be defined approximately, for finite A'^, 
as the point at which (r)a changes sign. This point can 
be easily read off, even for the small system sizes stud- 
ied here for which there is no rigorous relation with a 
true phase transition. As mentioned in the Introduction, 
the infinite-system value of Cc marks the transition to a 
spin glass phase. For the pSM in the thermodynamic 
limit, an implicit analytic expression for the critical en- 
ergy is known [27], yielding the approximate numerical 
value e-c ~ —1.17. The finite-system critical energies are 
plotted vs. the inverse system size 1/N in Fig. 3. De- 
spite the small system sizes, the results indicate a trend 
which appears to be nicely consistent with the analytic 
thermodynamic limit value. It is important to note that 
the finite N critical energies previously defined do not 
coincide with the ground state energy of the clusters. 



C. Index-resolved number of stationary points 

The index / of a stationary point a'^ , i.e., the number of 
negative eigenvalues of the Hesse matrix T-Lia^)^ is known 
to be a relevant quantity for the physical properties of the 
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FIG. 3. (color online) Finite-system critical energies ec, de- 
termined as the energy values at which the annealed average 
(r)a in Fig. 2 changes sign. For comparison, the thermody- 
namic limit value ec ~ —1.17 is plotted in red. 



system under consideration. In particular, minima are 
stationary points with / = and correspond to regions 
in phase space in vi^hich a trajectory may be trapped, 
which in turn is relevant for glassy dynamics. Higher- 
index stationary points, in contrast, allow the system to 
diffuse more easily due to their inherent instability. 

In Fig. 4 (top) we plotted, for several system sizes N, 
the (logarithm of) total number of stationary points vs. 
the index /. The number of stationary points has a max- 
imum at and is symmetric around / = N /2. This sym- 
metry is a direct consequence of the antisymmetry of 
the Hamiltonian (5) of the 3SM: For every stationary 



point & 



0' " 



0' "1' ■ 
rf,.. 



, fT^) with index /(cr*'), the point 
Pfj is also stationary and its in- 



— fT 

dex is /(P) ^N~I{cr 

For the pSM in the thermodynamic limit, analytic re- 
sults are available for the mean total number Afi of sta- 
tionary points of a given index /, 



lim —\nE( jVi) ^-ln2-- 



0.013, 



(17) 



see Eq. (2.21) of [27]. Interestingly, this result is indepen- 
dent of the index /. The numerical values for A// from 
our finite-system data are much larger than the infinite- 
system results, but the scaling behavior with N is roughly 
consistent with the expectation for small indices /. For 
/ > N/2, the curves in Fig. 4 (top) bend downwards as 
a consequence of the earlier mentioned symmetry of the 
3SM Hamiltonian (5). This strong finite-size effect can 
be reduced by considering the mean total number Afi of 
stationary points of a given index density i [see Fig. 4 
(bottom)]. The trend of the data in Fig. 4 (bottom), in- 
dicates a decrease of Afi with increasing system size, but 
our data do not allow an extrapolation to the thermody- 
namic limit. 

Similar to Fig. 4 (bottom), a unimodal shape of 
\nMi/N vs. i has previously been observed for several 
other models. Examples include the XY model with 
power law interactions for up to 13 degrees of free- 
dom [15]; the one-dimensional XY model, for which all 
the stationary points were found analytically for any N 
[28, 29]; the random-phase XY model in two dimensions. 
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FIG. 4. (color online) The total number of stationary points 
vs. index / (top) and vs. index density i (bottom) for various 
system sizes A'^. 



for which all stationary points were computed for small 
square lattices (up to a 3x3) [14, 22]; Lennard- Jones clus- 
ters for up to 14 atoms [30]. In Ref. [31], for large bulk 
systems consisting of weakly interacting sub-systems, Afi 
was analytically shown to follow a binomial distribution 
in /. 



D. Index-resolved complexity 

In Fig. 5 (top) the cumulative complexities Tj of sta- 
tionary points of index / = 0, 1, 3, and 10 are shown as 
functions of the energy density e for various system sizes. 
Not unexpectedly, lower-index stationary points are more 
numerous at low energies, higher-index stationary points 
at high energies. Below a certain threshold value around 
e = —0.9, only minima (stationary points of index / = 0) 
are present. 

Also for the index-resolved complexities F/ analytic 
results are known in the thermodynamic limit, see e.g. 
Eq. (2.16) of [27]. For a convenient comparison we show 
a plot of these results in Fig. 5 (bottom). For energies 
e > eth = — 2/-\/3 « —1.155 above the threshold energy, 
the cumulative complexities F/ in the thermodynamic 
limit are constant in e and independent of /, 



r/(e) = In 2/2 - 1/3 for all / and e > eth- 



(18) 



Below eth the inequality r/_i(e) > r/(e) holds for all 
/ = 1, . . . , A^. This inequality implies that, in the ther- 
modynamic limit, the complexity ro(e) of minima dom- 
inates the total complexity r(e) for all e < eth- Due 
to this dominance of minima, true ergodicity breaking 
occurs below this threshold, and the diffusivity of the 
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systems goes to zero. This makes the threshold energy 
eth a relevant characteristic of the dynamical behavior of 
the pSM [2, 4, 32]. 

Strictly speaking, the threshold energy eth is defined 
only in the thermodynamic limit where the heights of 
the barriers between stationary points diverge. For fi- 
nite N the barrier height remains finite, ergodicity is 
restored, and therefore one cannot expect a true diffu- 
sion arrest in the dynamics to take place at a particular 
energy. Nonetheless it would be interesting to study a 
finite-system analog of the threshold energy level in the 
pSM. For some models which in certain respects behave 
similar to the pSM, such a finite-system analog of eth can 
be defined as the energy below which a crossover takes 
place towards extremely slow dynamics. This slowing 
down is attributed to the suppression of escape direc- 
tions in the energy landscape [9, 10]. For the rather small 
system sizes at our disposal, we expect the relaxational 
dynamics of the system to be ruled by a spectrum of re- 
laxation time scales that are related to the iV-dependent 
heights of energy barriers in the low-energy regime dom- 
inated by minima. It would be interesting to analyze the 
dynamics of small clusters of spins in this regime and try 
to interpret the results in the light of the energy land- 
scape properties described here. This could shed light on 
the mechanism that leads, in the large-system limit, to 
the characteristic two-step relaxation of glassy systems. 

For the small system sizes we have studied, the finite- 
results do not even qualitatively resemble the thermody- 
namic limit results. For small indices / ^ iV, we at least 
observe the correct trend, as the large-e value of F/ is 
decreasing with increasing system size N. The trend is 
reversed (and inconsistent with the large-TV result) for 
larger values of /, but this finite-size behavior can be at- 
tributed to similar reasons as discussed at the end of Sec. 
IV C. 



E. Cumulative mean index density 




0.02 
0.01 



-0.01 
-0.02 




FIG. 5. (color online) Top: The cumulative complexities F/ 
of stationary points of index / = 0, 1, 3, and 10 as functions 
of the energy density e for system sizes A'^ — 14, 16, 18, and 
20. Among the curves of a given index 7, the energy ei at 
which Ti changes sign decreases with increasing system size, 
and this regularity allows one to easily relate system sizes to 
the various curves in the plot. Bottom: The corresponding 
thermodynamic limit results. 




One indication of the occurrence of a mode coupling 
transition at some value of the energy density e is, as 
mentioned in the Introduction, the vanishing of the mean 
index density i(e) [9, 10, 30, 33], as obtained by averag- 
ing the index density i = I/N of the Hesse matrix over 
all stationary points with energy e. For the reasons ex- 
plained above, we prefer to use a cumulative version of 
this quantity, 

{a": H{a^)<Ne} 

defined as the index density averaged over all stationary 
points with energy densities not larger than e. Results 
for J, averaged over 10 disorder realizations are shown for 
system sizes N = 14, . . . , 20 in Fig. 6. The cumulative 
mean index density is found to be a monotonically in- 
creasing function of the energy density e. saturating at a 



FIG. 6. (color online) Cumulative mean index density j vs. 
energy density e for various system sizes A'. 



value of 1/2 with increasing energy. This limiting value 
is again a consequence of the symmetries of the 3SM. 
Within the narrow range of system sizes we were able to 
probe, the A'^-dependence is very weak. 

As previously remarked, in the large-A'^ picture, the 
threshold energy eth below which the landscape is dom- 
inated by minima plays an important role. Naively one 
might try to identify the finite-system counterpart of eth 
with the value of e at which j{e) starts deviating from 
zero. Here we follow a different approach: In the large- 
A'^ limit and asymptotically for e > eth, the mean index 
density i{e) is expected to grow like a power law of the 
distance from the threshold energy (Eq. (55) of [33]). Ac- 
cordingly, a similar power law holds for the cumulative 
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FIG. 7. (color online) Finite-system threshold energies eth as 
obtained from fits to the cumulative mean index densities j(e) 
in Fig. 6. An extrapolation to small values of 1/N appears 
at least not inconsistent with the thermodynamic limit value 
eth = -2/^/3« -1.155. 



mean index density, 

j(e) (X (e- 



eth)" for e > eth 



(20) 



with exponent a ~ 3/2. Fitting such a power law to 
the data of Fig. 6, we find that the optimal fit param- 
eters a « 1.3 differs somewhat from the expected ther- 
modynamic limit value. The finite-system threshold en- 
ergies eth, though not an exact match, are at least in 
reasonable agreement with the large-A^ predictions. In 
Fig. 7, eth is plotted vs. inverse system size, and we ob- 
serve a correct trend towards the infinite-system value 
eth - -2/V3« -1.155. 



F. Hesse determinant 

In order to further characterize the stationary points of 
the Hamiltonian (5), we computed the determinant of the 
Hesse matrix at each of the stationary points. Similar, 
but complementary, to the index of a stationary point, 
also the Hesse determinant represents a way of projecting 
the information contained in the Hesse matrix onto a 
single number. In energy landscape studies, the Hesse 
determinant has been used previously in the study of 
phase transitions, also in the absence of disorder. One 
of the key results in this context relates the vanishing of 
the Hesse determinant in the thermodynamic limit to the 
occurrence of a phase transition [34, 35]. More precisely, 
under some technical assumption, in these references it 
was shown that, if there exists a sequence of stationary 
points (T/„^ for iV G N such that the energy density of 



'{N) 

this sequence converges. 



N 



hm ff(a^^))/iV, 



and such that 



N 



lim |detH(cr^^J|^^^ = 



(21) 



(22) 



then the stationary points may induce a phase transition 
at the energy density ec in the thermodynamic limit. Eq. 
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FIG. 8. (color online) For all stationary points a" of our 10 
disorder realizations of the 3SM, the pairs (e, D) are plot- 
ted, where e = H{a^)/N is the energy density and D = 

|det'H(o-°)|^''^ is the rescaled Hesse determinant at a station- 
ary point. The top plot is for system size A*' = 16, the bottom 
one for N = 20. 



22 implies that stationary points have to become suf- 
ficiently "flat" in the thermodynamic limit in order to 
be capable of inducing a phase transition. This scenario 
was confirmed in a number of exactly solvable long-range 
interacting models (mostly with mean-field interactions) 
[34-36] , and it also led to the analytic prediction of the 
exact transition energy of the self-gravitating ring model 
[37]. 

On the basis of finite-system data, it is of course 
not possible to construct infinite sequences of stationary 
points (T^^j for N £ N. But we can try and see whether, 
for the small system sizes at our disposal, the Hesse de- 
terminant exhibits a tendency to vanish at some value 
of the energy, possibly close to the threshold energy or 
critical energy of the 3SM. To investigate this, we have 
computed Hesse determinants and energy densities for 
all stationary points of all our 10 disorder realizations of 
the model. The data, plotted in Fig. 8, do not reveal any 
indication of a vanishing Hesse determinant. Moreover, 
comparing the data for system sizes iV = 16 (top) and 

= 20 (bottom) in Fig. 8, no significant size-dependence 
is visible. 

This is a negative result, and its interpretation is not 
entirely clear. On the one hand, the absence of any sig- 
nificant differences between the plots for A^ = 16 and 
A^ = 20 may be seen as an indication that the range of 
system sizes we are able to deal with is way too small to 
observe finite-size scaling effects that allow for an ex- 
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trapolation towards the thermodynamic hmit. But it 
may of course also be that, even for much larger systems, 
the Hesse determinant remains bounded away from zero. 
Note that such a scenario would not be in contradiction 
with the criterion (22): This criterion specifies necessary 
conditions for a sequence of stationary points to induce 
a phase transition. And although examples are known 
where such stationary points indeed are responsible for 
the occurrence of a phase transition, this need not be the 
case for other models. And indeed, several model systems 
are known in which stationary points are not related to 
the occurrence of a phase transition [17, 38-44], and this 
might also be the case for the 3SM. 

In summary, the exponentially (in N) large number of 
stationary points renders the criterion (22) very imprac- 
tical for a naive numerical approach like ours. Only with 
intuition about the relevant types of stationary points, 
and a numerical algorithm focusing exclusively on those 
points, might a numerical test of the condition (22) be 
viable. 



V. DISCUSSION AND CONCLUSIONS 

Wc have studied the potential energy landscape of the 
mean-field p-spin spherical model (pSM) for p = 3. By 
means of the numerical polynomial homotopy continu- 
ation method, we obtained a complete characterization 
of the energy landscape, in the sense of computing all 
stationary points of the Hamiltonian (5) for a number 
of disorder realizations. Since the number of stationary 
points grows exponentially with the system size the 
capacity to compute all the stationary points is simulta- 
neously a virtue and a limiting factor, as it restricts the 
method to rather small system sizes (up to = 20 in 
our study). 

The pSM is known as a prototype of a fragile glass 
forming system and for several of its landscape character- 
istics, like the number of stationary points, indices, and 
associated complexities, analytical results are known in 
the thermodynamic limit. The research reported in the 
present article complements these thermodynamic limit 
results with their finite-system counterparts. The aim of 
such a study is to gain insight into how certain infinite- 
system properties of the model emerge from finite N clus- 
ters. In particular, when dealing with mean field models 
which show artificially divergent energy barriers when 
N ^ oo, the finite-size scaling of the barrier height with 
N is a key quantity for understanding the slow dynam- 
ics of the model. For the small system sizes studied, 
we found results that, though roughly compatible with 
large- A^ results, show significant finite-size effects. 

The index- resolved number of stationary points. A//, 
shows a maximum at / = N/2 and, due to a symmetry 
of the Hamiltonian (5), is symmetric around this point. 
At variance with this result, analytical calculations in the 
thermodynamic limit give just a constant value, indepen- 
dent of /. While the change with increasing A^ is roughly 



compatible with the asymptotic result, we can not draw a 
strong conclusion about this point which deserves further 
study with larger systems. 

When studying the cumulative complexities r(e), we 
found it convenient to focus on the annealed disorder av- 
erages as defined in (16). These quantities have the ad- 
vantage of changing sign at a certain value of the energy 
density e, which we define as the finite-A 'pseudo'-critical 
energies. These can be seen as finite-system counterparts 
of the critical energy in the thermodynamic limit. Al- 
though, for the small sizes studied, the actual values of 
the pseudo-critical energies are not close to the thermo- 
dynamic limit results, the trend is the correct one. An 
interesting outcome of these results is that the number 
of saddles is exponential in the system size for energies 
down to a critical value (N dependent) which is always 
larger than the ground state energy, i.e., this important 
property is not limited to the thermodynamic limit. 

Strong finite size effects are also seen in the index- 
resolved cumulative complexities. Nevertheless, in agree- 
ment with the largc-A picture, the low energy sector is 
found to be dominated by small-index stationary points. 
In particular, at the lowest energies, only stationary 
points of index / = (i.e., minima) occur. Upon in- 
creasing the energy, stationary points with higher in- 
dices appear and gradually become more numerous than 
lower-index stationary points. This should also reflect in 
the relaxation properties of the model. In the large-A' 
limit, minima dominate completely for energy densities 
between the critical one Cc and the threshold energy eth- 
This leads to the well known two-step relaxation of the 
model in that energy regime. For small clusters we ex- 
pect that a hierarchy of time scales will be present, which 
reflects the stability properties of the dominant set of sta- 
tionary points at each energy. 

Different from the quantities mentioned earlier, the 
mean index density shows rather weak flnite-size effects. 
Its behavior is qualitatively similar to the large A^ results, 
rising like a power law above the 'pseudo'-threshold en- 
ergy. From fits to our data we obtain an exponent a w 1 .3 
for such a power law, reasonably close to the asymptotic 
result a = 3/2. From the same fits finite-A^ approxima- 
tions to the threshold energy can be extracted. These 
energies show the correct trend with increasing A^ and 
are systematically larger than the corresponding pseudo- 
critical ones. The energy intervals so defined should be 
the small-A^ analogs of the 'minima-dominated' regime 
known to be present in the thermodynamic limit. If this 
has implications to the finite N dynamics is not known 
and is a point that certainly deserves to be studied. 

Finally wc studied the rescaled Hesse determinant D, 
evaluated at the stationary points of the 3SM. It is 
known, at least for certain models studied in the liter- 
ature, that the vanishing of this quantity at some value 
of e in the thermodynamic limit may signal the pres- 
ence of a thermodynamic phase transition at that energy 
density. Despite the fact that the model shows a phase 
transition at Cc in the thermodynamic limit, wc found no 
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evidence of a tendency to zero of the Hesse determinant 
at any value of e. There are many possible reasons for 
this to happen, the small system-sizes considered being 
one of them. Another reason may lie in the determinant 
criterion itself: It is based on the assumption that sta- 
tionary points are 'at the origin' of a phase transition. 
Although this hypothesis is correct in many models, it is 
known to fail in others [17, 38-44], and in this case the 
reasoning behind the criterion breaks down. 

The present paper focuses on properties of the poten- 
tial energy landscape only. Since many of these proper- 
ties are known to reflect also in dynamical and thermo- 
dynamical properties of finite systems, a study of these 
latter properties and a comparison to finite-TV energy 
landscape properties may lead to further insights on the 
interplay between landscape and dynamics. This is a 



project worth being addressed in the future. 



ACKNOWLEDGMENTS 

D.M. was supported by the U.S. Department of En- 
ergy under Contract No. DE-FG02-85ER40237, and he 
would like to thank Tarek Anous and Frederik Denef for 
introducing him to the p-spin model. The numerics were 
carried out on Fermilab USQCD cluster. M.K. acknowl- 
edges support by the Incentive Funding for Rated Re- 
searchers programme of the National Research Founda- 
tion of South Africa. D.A.S. acknowledges partial sup- 
port from a research fellowship by CNPq, Brazil. 



[1] A. Crisanti and H.-J. Sommers, Z. Phys. B 87, 341 
(1992). 

[2] J. Kurchan, G. Parisi, and M. A. Virasoro, J. Phys. I 

France 3, 1819 (1993). 
[3] T. Castellani and A. Cavagna, J. Stat. Mech. 2005, 

P05012 (2005). 

[4] A. Crisanti and H.-J. Sommers, J. Phys. I France 5, 805 
(1995). 

[5] A. Cavagna, I. Giardina, and G. Parisi, Phys. Rev. B 

57, 11251 (1998). 
[6] A. Crisanti, L. Leuzzi, and T. Rizzo, Eur. Phys. J. B 36, 

129 (2003). 

[7] A. Cavagna, Europhys. Lett. 53, 490 (2001). 

[8] J. P. Bouchaud, L. Cugliandolo, J. Kurchan, and 

M. Mezard, Physica A 226, 243 (1996). 
[9] T. S. Grigera, A. Cavagna, I. Giardina, and G. Parisi, 
Phys. Rev. Lett. 88, 055502 (2002). 
[10] L. Angelani, R. Di Leonardo, G. Ruocco, A. Scala, and 

F. Sciortino, Phys. Rev. Lett. 85, 5356 (2000). 
[11] P. G. Debenedetti and F. H. Stillinger, Nature 410, 259 
(2001). 

[12] A. J. Sommese and C. W. Wampler, The Numerical So- 
lution of Systems of Polynomials Arising in Engineering 
And Science (World Scientific, Singapore, 2005). 

[13] D. Mehta, Lattice vs. Continuum: Landau Gauge Fixing 
and 't Hooft-Polyakov Monopoles, Ph.D. thesis. The Uni- 
versity of Adelaide (2009), Australasian Digital Theses 
Program. 

[14] D. Mehta, A. Sternbeck, L. von Smekal, and A. G. 

WiUiams, PoS QCD-TNT09, 25 (2009). 
[15] D. Mehta, Phys. Rev. E 84, 025702(R) (2011). 
[16] D. Mehta, Adv. High Energy Pliys. 2011, 263937 (2011). 
[17] M. Kastner and D. Mehta, Pliys. Rev. Lett. 107, 160602 

(2011). 

[18] R. Ncrattini, M. Kastner, D. Mehta, and L. Casetti, 
"Exploring the energy landscape of XY models," 
arXiv:1211.4800. 

[19] D. Mehta, J. D. Hauenstein, and M. Kastner, Phys. Rev. 
E 85, 061103 (2012). 

[20] M. Maniatis and D. Mehta, "Minimizing Higgs poten- 
tials via numerical polynomial homotopy continuation," 
arXiv:1203.0409. 



[21] D. Mehta, Y.-H. He, and J. D. Hauenstein, JHEP 1207, 
18 (2012). 

[22] C. Hughes, D. Mehta, and J.-I. SkuUerud, Ann. Physics 
331, 188 (2013). 

[23] D. Martinez-Pedrera, D. Mehta, M. Rummel, and 
A. Westphal, "Finding all flux vacua in an explicit ex- 
ample," arXiv: 1212.4530. 

[24] Y.-H. He, D. Mehta, M. Niemerg, M. Rummel, and 
A. Valeanu, "Exploring the potential energy landscape 
over a large parameter-space," arXiv:1301.0946. 

[25] J. D. Hauenstein and F. Sottile, ACM T. Math. Software 
38, 28 (2012). 

[26] D. Mehta, J. D. Hauenstein, and D. J. Wales, (2013), 

arXiv: 1302.6265 [cond-mat .stat-mech] . 
[27] A. Auflinger, G. Ben Arous, and J. Cerny, Commun. 

Pure Appl. Math. 66, 165 (2011). 
[28] D. Mehta and M. Kastner, Annals Phys. 326, 1425 

(2011). 

[29] L. von Smekal, D. Mehta, A. Sternbeck, and A. G. 

Williams, PoS LAT2007, 382 (2007). 
[30] J. P. K. Doye and D. J. Wales, J. Chem. Phys. 116, 3777 

(2002) . 

[31] D. J. Wales and J. P. K. Doye, J. Chem. Phys. 119, 12409 

(2003) . 

[32] L. F. Cugliandolo and J. Kurchan, Phys. Rev. Lett. 71, 
173 (1993). 

[33] J. Kurchan and L. Laloux, J. Phys. A: Math. Gen. 29, 
1929 (1996). 

[34] M. Kastner and O. Schnetz, Phys. Rev. Lett. 100, 160601 

(2008) . 

[35] M. Kastner, O. Schnetz, and S. Schreiber, J. Stat. Mech. 

2008, P04025 (2008). 
[36] M. Kastner, Phys. Rev. E 83, 031114 (2011). 
[37] C. Nardini and L. Casetti, Phys. Rev. E 80, 060103 

(2009) . 

[38] F. Baroni, Transizioni di fase e topologia dello spazio del- 
le configuraztom di modelli di campo medio, Master's the- 
sis, Universita degli Studi di Firenze (2002). 

[39] D. A. Garanin, R. Schilling, and A. Scala, Phys. Rev. E 
70, 036125 (2004). 

[40] A. Andronico, L. Angelani, G. Ruocco, and F. Zamponi, 
Phys. Rev. E 70, 041101 (2004). 



10 



[41] M. Kastner, Phys. Rev. Lett. 93, 150601 (2004). [44] L. Angelani and G. Ruocco, Phys. Rev. E 76, 051119 

[42] L. Angelani, G. Ruocco, and F. Zamponi, Phys. Rev. E (2007). 

72, 016122 (2005). 
[43] I. Hahn and M. Kastner, Phys. Rev. E 72, 056134 (2005). 



